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ABSTRACT 

We present an exploration of the significance of Carbon/Oxygen phase separation in white dwarf stars in the 
context of self-consistent evolutionary calculations. Because phase separation can potentially increase the calcu- 
lated ages of the oldest white dwarfs, it can affect the age of the Galactic disk as derived from the downturn in the 
white dwarf luminosity function. We find that the largest possible increase in ages due to phase separation is ~ 1.5 
Gyr, with a most likely value of approximately 0.6 Gyr, depending on the parameters of our white dwarf models. 

The most important factors influencing the size of this delay are the total stellar mass, the initial composition 
profile, and the phase diagram assumed for crystallization. We find a maximum age delay in models with masses 
of ~ O.6M , which is near the peak in the observed white dwarf mass distribution. In addition, we note that 
the prescription which we have adopted for the mixing during crystallization provides an upper bound for the 
efficiency of this process, and hence a maximum for the age delays. More realistic treatments of the mixing 
process may reduce the size of the effect. We find that varying the opacities (via the metallicity) has little effect 
on the calculated age delays. 

In the context of Galactic evolution, age estimates for the oldest Galactic globular clusters range from 11.5 to 
16 Gyr, and depend on a variety of parameters. In addition, a 4 to 6 Gyr delay is expected between the formation 
of the globular clusters and that of the Galactic thin disk, while the observed white dwarf luminosity function 
gives an age estimate for the thin disk of 9.5^ g Gyr, without including the effect of phase separation. Using the 
above numbers, we see that phase separation could add between to 3 Gyr to the white dwarf ages and still be 
consistent with the overall picture of Galaxy formation. Our calculated maximum value of <; 1.5 Gyr fits within 
these bounds, as does our best guess value of ~ 0.6 Gyr. 

Subject headings: dense matter — equation of state — stars: evolution — white dwarfs 



1. ASTROPHYSICAL CONTEXT 

The phenomenon of phase separation and crystallization ex- 
ists within the larger context of whit e dwarf cooling. Since the 
time of Mestel's original treatment ( Mestel 1952), much work 
has been done, both to improve the input physics of the models 
and to make more complete observat ions of the wh ite dwarf lu- 
minosity function (WDLF). In 1987, Winget et al. showed that 
the observed downturn in the WDLF could be understood in 
terms of a finite age for the Galactic disk, and that the WDLF 
could therefore in principle be used to determine an age for 
the local Galactic disk. Using the preliminary results from 



the Oswalt et al 



(1996) sample gives an age of 9.5^q§ 



Gyr. 



Taki ng the error estim ates at face value, these results differ by 
2a. Wood & Oswalt (1998) conducted Monte Carlo simula- 



tions and found that it is unlikely that both samples are consis- 
tent with the same parent population. Further investigation will 
be needed to resolve the cause of this discrepancy. 

In addition to the uncertainties in the observed WDLF, the 
way we treat various physical processes in white dwarf interi- 
ors greatly affects the ages that we derive for them. After the 
prediction in the early 1960's that white dwarf s should un dergo 
a phase tr ansitio n and crys tallize as they cool ( A.brikosov 1960; 



Liebert, Dahn, & Monet (1988, hereafter LDM) for the ob- 



served WDLF, they obtained an age for the local Galactic disk 
in the range 7-10 Gyr. Since then, Wood has made more de- 
tailed calculations using improved input physics, Galactic evo- 
lution models , and WD parameters to constrain this age even 
further |Wood| (1990, 1992, 19 95). Historicall y, these develop 



Kirzhnits 1960; Salpetei 1961), Mestel & Ruderman (1967) 



and Van Horn (1968) estimated that the associated release of 



latent heat during this process would be la rge enough to delay 



ments were foreshadowed by Schwarzschild (1958), Schmidl 
(1959), and p'Antona & Mazzitelli| (1978), all of whom con- 
sidered white dwarf evolution in a Galactic context. 

Two observational surveys within the last ten years stand out 
in their importance to the field. First, Liebert et al. (1988) pro- 
duced a WDLF containing 43 cool field WD's, which was the 
largest such sample size up to that point in time. More recently, 
Oswalt et al.| (1996) produced a WDL F of 50 cool WD's in wide 
binaries. Using the models of Wood the LD M sam ple yields an 
age for the Galactic Disk of ~ 7.5 ± 1 Gyr (|Wood| 1995), while 



the cooling of white dwarfs significantly. Lamb & Van Horn 
(1975) included this energy release as part of their evolutionary 
calculations o f a 1 M Q pure carbon white dwarf. 

SteVenUon (1977) was the first to propose a phase separation 
model that might affect white dwarf cooling times by provid- 
ing an additional source of energy analogous to the release of 
latent heat. This model ha d a carbon core with trace amounts 
of iron. In a later model, Stevenson (1980) suggested that a 
uniform mixture of carbon and oxygen would become chemi- 
cally differentiated as a result of the crystallization process. Be- 
cause such a redistribution of elements could lower the binding 
(non-thermal) energy of the star, the change in energy would be 
added to the thermal energy, and hence the luminosity, of the 
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star. This would increase the time for a white dwarf to cool to 
a given luminosity, and would extend the apparent age of the 
Galactic disk as derived from the WDLF. 

Estimates of the amount by which the age of the local Galac- 
t ic disk might be exte nded have ranged from 0.5 Gyr to 6 Gyr 
(|Mochkovitch| 1983 ; |Barrat, Hansen, & Mochkovitch] 1988; 
Garcfa-Berro et al.| 1988; IChabrier et al.| 1993; Se gretain & 
Chabrier 1993; |Hernanz et alj 1994; gegretain et al.| 1994; I s- 
ern et al. 1997; [Salaris et alj 1997), although rec ent estimates 



pressure, this C-e nhanced fluid l ayer i s mixed via a Rayleigh 



Taylor instability (Mochkovitch 1983; Isern et al. 1997) with 



have been on the smaller end of this range, e.g., Salaris et al 
(1997) calculate a delay of ~ 1.0 Gyr. Most of this spread in 
calculated age delays comes from differences in the assumed 
phase diagram, although the assumed C/O profile also has a 
large effect. 

In the context of Galactic evolution, age estimates for the 
olde st Galactic globular clusters range from 13-16 Gyr (P ont 
et al. 1998) to 11.5 ~ 1.3 Gyr ( piaboyeret alj 1998), and de- 
pend on a variety of parameters. In addition, a 4 to 6 Gyr delay 
is expected between the formation of the globular clusters and 
that of the Galactic thin disk (e.g., Burkert, Truran, & Henslei 



1992; phiappini, Matteucci, & Grattorj 1997), while the ob- 
served white dwarf luminosity f unction gives an age estimate 

1996), without 



for the thin disk of 9.5^ Gyr ( pswalt et al 



including the effect of phase separation. Using the above num- 
bers, we see that phase separation could add anywhere from 
to 3 Gyr to the white dwarf ages and still be consistent with the 
overall picture of Galaxy formation. 

In this paper, we examine the sensitivity of this calculated 
age delay to the various physical assumptions by varying the 
initial C/O profile of the white dwarf models, their total mass, 
and their H and He surface layer masses. In addition, we exam- 
ine the effect of using two different pub lished phase diagrams 
for the phase sepa ration process, that of Segretain & Chabriei 
(1993) and that of |lchimaru, Iyetomi, & Ogata| (1988). 

Our work improves upon previous calculations of the age de- 
lay in that we use self-consistent evolutionary models. In par- 
ticular, ^urjno^ds_use_tfiejno^ernOr^L^ & 
Rogers 1993) instead of the older Cox-Stewart opacities, and 
we are able to treat self-consistently the age delay as a function 
of total stellar mass, instead of using a relation scaled by mass 
for the connection between the core temperature and the sur- 
face luminosity. Finally, we are able to examine surface layer 
masses suggested by more recent asteroseismological investi- 
gations ( ]Clemens| 1993). 



2. THE PHYSICS OF PHASE SEPARATION 
2.1. Chemical Redistribution 

Our present physical picture for the phenomenon of phase 
separation in white dwarf stars is as follows. As a white dwarf 
cools, it eventually reaches a temperature when its central re- 
gions begin to crystallize. This occurs when the thermal en- 
ergy of the ions becomes much smaller than the energy of the 
Coulomb interactions between neighboring ions. As a result, 
the ions settle into lattice sites and lose the ability to move 
freely in three dimensions. 

If the white dwarf interior is initially a mixture of C and O, 
then recent calculations indicate that the solid which crystal- 
lizes will have a higher O conte nt than the fluid from which 
it formed (Ichimaru et al. 1988; Segretain & Chabriei 1993). 
Thus the crystallizing region of the white dwarf becomes O- 
enhanced and the fluid layer overlying this region becomes C- 
enhanced. Since the C is slightly less dense than the O at a given 



the layers above, and C is transported outward from the cen- 
ter. As the white dwarf continues to crystallize, the O-enhanced 
crystalline core also continues to grow, with the net result that 
O is transported inward in the white dwarf and C is transported 
outward. Thus, the chemical composition profile after signif- 
icant crystallization has occurred is different from the profile 
before crystallization. 

Just how different this profile is depends on the particular 
phase diagram which is adopted for the process. In a "spin- 
dle" diagram, the solid which forms always has an enhanced 
concentration of the higher charge element (in this case oxy- 
gen), and the temperature of crystallization of the mixture lies 
between that of the individual elements. An "azeotropic" dia- 
gram differs from this in that there is a range of concentrations 
for which crystallization takes place below the temperature of 
crystallization of either of the pure elements. This is some- 
what analogous to the phenomenon of "supercooling." Finally, 
a "eutectic" phase diagram is one in which there is a near total 
separation of the higher and lower charged ions upon crystal- 
lization, resulting in a segregation of th e two chem ical species. 

Stevenson's original phase diagram (Stevenson 1980) was a 
eutectic phase diagram with C and O being immiscible in the 
solid phase, with the result that a pure O core would be formed 
in the models during crystallization. Using a density functional 
approach, Barrat et al. (1988) calculated a phase diagram of 
spindle type. In this case, the solid which forms is a C/O alloy, 
but with the O content of the solid enhanced relative to that of 
the fluid out of which it formed. 

This problem was revisited by Ichimaru et al. (1988). They 
found that Stevenson's initial prediction of a eutectic phase 
diagram was an artifact of his use of the random-alloy mix- 
ing (RAM) model for the internal energies in the solid phase. 
By comparison with Monte Carlo simulations, they found that 
the linear mixing formula is more accurate for the solid phase. 
They then used density-functional theory to derive a phase di- 
agram of azeotropic type, which is shown as the dashed line in 
Figure [|. This diagram is similar to the spindle diagram, with 
the exception that there is a range of compositions for which 
the crystallization temperature is less than the crystallization 
temperature for e ither of the pure compo sitions. 

Most recently, Segretain & Chabrier (1993) used a density- 
functional approach to derive phase diagrams for arbitrary 
binary-ionic mixtures, as a function of Z1/Z2, where Z\ and 
Z2 are the nuclear charges of the two chemical species. For C 
and O (Z1/Z2 = 0.75), they obtain a phase diagram of spindle 
type, which is shown as the solid line in Figure]]]. 

As shown in Figure [lj, these diagrams of Ichimaru et al. and 



Segretain & Chabriei] differ slightly in the composition changes 



during crystallization, as well as in the temperatures at which 
crystallization takes place. As a result, they produce different 
chemical composition profiles after crystallization and different 
age delays. 

We mention that the validity of the azeotropic phase diagram 
of Ichimaru et al. (1988) has rec ently been called into question 
by |DeWitt, Slattery, & Chabrier| (1996). DeWitt et al. find that 
the azeotropic diagram o f Ichimaru et al. res ults from their use 
of the fitting function of |Ogata & Ichimaru (1987) for the en- 
ergies of the pure phases. The slight inaccuracies in this fitting 
function lead to incorrect values of the free energy for the pure 
phase, which in turn leads to a spurious departure from the lin- 
ear mixing rule in the liquid. With this proviso, we have elected 
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FIG. 1 . — Phase diagrams for a C/O mixture as computed by Ichimaru et al. (1988, dashed line) and Segretain & Chabriei (1993, solid line), where the vertical 
axis is in units of the crystallization temperature of C, 7c, and the horizontal axis is the C mass-fraction. The solid line is of "spindle" type, while the dashed line is 
that of an "azeotrope." The principal feature of the azeotrope is that there is a range of compositions for which the crystallization temperature of the mixture is less 
than that of either of the two constituents in the pure state. 



to include calculations based on the phase diagram of Ichimaru 
et al. for purposes of comparison, in order to illustrate some of 
the uncertainties to which these calculations can be subject. 

2.2. Energy Release 

Because the distribution of C and O within a model changes 
during phase separation, the density profile changes as well. At 
a given pressure, O is slightly denser than C. A model which 
has undergone phase separation has more oxygen in its core, 
and thus a slightly larger concentration of mass in its central 
regions. As a result, the phase-separated model is more tightly 
bound gravitationally. 

While it may be convenient to think of the energy which is 
released as being due solely to the change in the gravitational 
potential energy of the star, this is only part of the story. The 
relevant quantity is actually the total binding energy of the star, 
£bind, which is the sum of all the nonthermal (structural) sources 
of energy. As such, it acts as a potential energy for the configu- 
ration. £bind can b e written as 



^bind — ^grav ^de° ^coul 



(1) 



where £grav,£deg, and iScoui are the respective energy contribu- 
tions from gravitational interactions, kinetic energies of the de- 
generate electrons, and Coulomb interactions among the differ- 
ent charged particles (ions and electrons). 

As phase separation occurs and the central regions become 
oxygen enriched, the central density of the model increases. 
Thus, £g rav becomes more negative, as does E cou \. Ed eg , how- 
ever, becomes more positive, since the Fermi energy of the elec- 
trons increases with increasing density. Summing these contri- 
butions, we find that there is a net decrease in Zibind f° r the mod- 
els considered here. Due to conservation of energy, this energy 



must be used to increase the thermal energy of the ions, which 
are the only significant repository of thermal energy in the cores 
of our white dwarf models. This energy, then, is available to be 
radiated away and acts as an additional luminosity source. 

The various contributions to the photon luminosity L and the 
neutrino lumin osity L v of the white dwa rf may be formally 



written as (e.g., [sern et al. 1997; Hhabrier 1998) 



L+L v -- 



Mm 



dm 



dT 



dTjyy dt 



du 
dX 



T.V 



dXo 
dt 



(2) 



where V = l/p is the specific volume, Xo is the mass-fraction 
of the heavier of the two chemical species (in this case, oxy- 
gen), and u is the internal energy per unit mass, which contains 
thermal, electron degeneracy, and Coulomb contributions. The 
first term in the integrand on the righthand side of equation || is 
due to the heat capacity of the core, which includes the release 
of latent heat of crystallization, while the second term gives 
the contribution to the luminosity due to volume changes, and 
is usually small in white dwarfs since the pressure P is only 
a weak function of the temperature. The final term gives the 
luminosity due to the changing chemical composition profiles 
within the white dwarf. This is the term we will study in our 
numerical calculations. 

As a check on the direct evolutionary calculations, we can es- 
timate the age delay produced by a given energy release. If we 
denote by dE a small amount of energy which is released dur- 
ing the process of phase separation, and if we assume that this 
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energy is quickly radiated, then we can calculate an estimated 
age delay f^: 

f dE 

tj= / — • (3) 



In the context of a sequence of evolutionary models, this inte- 
gral is operationally a sum, since a given model is computed at 
discrete points in time, luminosity, etc. Furthermore, since the 
energy AEj is released between luminosities L,_i and L,-, say, 
the average luminosity at which the energy is released is ap- 
proximately (L,_i +L,)/2, so the discrete version of equation || 
becomes 

^Ettt^tt- (4) 



We have used equation^ as an alternate prescription to calculate 
age delays. For the larger energy releases, td computed in this 
way agrees with the delay calculated from the self-consistent 
evolutionary calculations, and for small energy releases it pro- 
vides a better estimate since the small energies can become 
masked in the numerical noise of the evolutionary calculations. 

3. NUMERICS 

The basis for these calculations is WDEC, the White Dwarf 



Evoluti onary C ode, as described in Lamb & Van Horn (1975), 
and in Wood (1990). Our current versi on use s the u pdated 
OPAL opacity tables ( jlglesias & Rogers] 1993; |Wood| 1993). 
We use the additive volume technique to treat the equation of 
state of the C/O mixture in the cores of our models. 

3.1. The Melting Curve 

Our criterion for crystallization is given by the phase dia- 
gram which we adopt, with the following caveat. Our equation 
of st ate (EOS) is based on the Lamb EOS code (L amb & Van 
Horn 1975), which has T ~ 160 at crystallization. Here, T = 
Z 2 e 2 1 (r)kBT is the ratio of Coulomb energy between neighbor- 
ing ions to each ion's k inetic energy. More recent calculations 
indicate that T ~ 180 ( |Ogata & Ichimaru| 1987). As a result, 
our values for the crystallization temperature of C, 7c. x tai, are 
too high by a factor of - 180/160= 1.125. 

To remedy this situation, we could simply adjust 7c. xta i down- 
ward accordingly, and we have done this for a few runs. This 
is inconvenient, however, because it places us at the edge of 
our EOS tables which were calculated with T ~ 160. Instead, 
we apply a correction factor to our calculated age delays which 
takes into account the fact that crystallization/phase separation 
occurs at lower central temperatures, and therefore lower lumi- 
nosities, than is calculated directly in our models. This correc- 
tion to the calculated age delays is typically of order 25%. For 
example, an age delay due to phase separation computed with 
r ~ 160 might be ~ 1 Gyr, which with the more physical value 
of r ~ 180 would be ~ 1.25 Gyr. We find that this procedure 
of computing the age delays based on the corrected luminosity 
during crystallization to be accurate to within 1-2%. 

3.2. Implementation in WDEC 

The calculation of the evolutionary sequences is "quasi- 
static" in the sense that we compute a sequence of static models 
separated by finite steps in time. Each static model represents 
the cooling white dwarf at a different age and luminosity. We 
include the physics of phase separation using the same approxi- 
mation: we assume that the timescale for any mixing which oc- 
curs is short compared to the individual evolutionary timesteps 



(see section 3__3 of this paper; Mochkovitcl 1983; Isern et al] 
1997), and we assume that the binding energy which is released 
by this process can be modeled by some su itably chosen local 
energy generation rate, e ps (e.g., Isern et al. 1997). 

The phase separation calculation may therefore be broken 
into three sections. The first part involves obtaining the chang- 
ing composition profile as a function of the crystallized mass- 
fraction, while the second part is the calculation of the cumula- 
tive energy released, also as a function of the crystallized mass- 
fraction. The final part is the calculation of the value of e ps , 
which is the energy locally deposited per unit mass per unit 
time. Our implementation of the complete problem is self- 
consistent in that we let e ps vary as the compositional profile 
changes due to crystallization, as WDEC iterates to a converged 
model. 

The first part of the overall problem relates to the compo- 
sition of the crystalli zing layers. Using the phase diagram of 
Segretain & Chabrier (1993) or Ichimaru et al. (1988), we com- 
pute the final composition profile of the model given the ini- 
tial profile, before doing a full evolutionary calculation. This 
is possible because the composition of the crystals which are 
forming is determined solely by the mass fractions of C and O 
which are present in the fluid phase, and not by the tempera- 
ture and density of the medium (the temperature and density of 
course determine when the fluid crystallizes, but given that it 
is crystallizing, the chemical composition of the solid is deter- 
mined solely by the composition of the fluid). We are therefore 
required to compute only once, at the onset of crystallization, 
the composition profile as a function of the crystallized mass 
fraction. At subsequent evolutionary times, we use this relation 
and the current crystallized mass-fraction to interpolate onto the 
composition grid, which is a computationally convenient proce- 
dure. 

We take this same approach for the calculation of the energy 
released. At the onset of crystallization, we calculate the total 
amount of energy released as a function of M xt - d \/M+, using the 
relation 

SE =L (ad./**' <5) 

where SE is the binding energy released by the composition 
change SXo- Since these changes in composition are with re- 
spect to the pre-crystallization state, we are in effect holding 
both the temperature and density profiles constant for all sub- 
sequent phases of crystallization. Holding the temperature pro- 
file constant is a quite reasonable approximation, since the vast 
majority of the mass in the the white dwarf model is strongly 
degenerate for the temperature range of interest. Similarly, we 
expect the changes in the density profile to be small (^ <J 1%) 
even in the presence of composition changes. This is a con- 
sequence of the fact that the equations of state for carbon and 
oxygen are very similar in the strongly degenerate regime, i.e., 
/i e , the atomic mass per electron, is 2.0 for both elements. This 
suggests that this approach would not necessarily be as accurate 
for carbon and iron, since /i e = 2.15 for iron. 

Figure ^ provides the final justification for our assumptions. 
The filled dots represent the energy released as computed self- 
consistently at each evolutionary time step, and the solid line 
is the calculated energy released assuming a static density and 
temperature profile as described in the above paragraph. The 
best agreement is for smaller amounts of crystallization, since 
these models differ the least from the initial static model. Even 
near complete crystallization, however, the difference between 
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the two values is less than 0.5%, justifying our assumptions. 
Computationally, it is very convenient to compute the energy 
release just once at the outset and then interpolate using the 
present value of the crystallized mass fraction. This allows 
WDEC to avoid doing a calculation of the energy release for 
each iteration of each model, which would significantly affect 
the speed of the calculations. 

Because all our calculations are done on evolutionary 
timescales, we do not have any information about the actual 
dissipative processes which are responsible for depositing the 
energy of phase separation locally. Indeed, without an accurate 
hydrodynamic model of the mixing process, this is not possi- 
ble. Fortunately, it is more important to know the total energy 
released rather than exactly how this energy is deposited within 
the white dwarf model. This is because the core has a very 
high thermal conductivity, which tends to smooth out the tem- 
perature distribution. Thus, wherever the energy is initially de- 
posited, it will soon be shared throughout the co re; inde ed, an 
isothermal core was an assumption of the original |Mestel theory 
(1952), and it is still a very accurate description of the p hysics 



in the interiors o f white dwarfs (e.g., Garcfa-Berro et al 
Segretain et al.| 1994) 



We therefore choose e ps such that the 
local temperature is increased by the same fractional amount 
throughout the core, i.e., 4j£ = const., while we simultaneously 
require that the total energy deposited in this way is equal to 
the energy released due to phase separation in a given timestep. 
Th is is somew hat analogous to the analytical approach outlined 
by Isern et al. (1997), although we developed our approach for 
ease of numerical implementation. 

There is one final adjustment which we make to the value of 
e ps as calculated above. It is due to the fact that WDEC calcu- 
lates models quasi-statically, so that e ps is assumed to have been 
constant during the last time step taken, when in fact it may 
have changed by a substantial amount. Put another way, the 
value of e ps which WDEC calculates should be associated with 
the average luminosity of the present and previous timesteps, 
not just the current luminosity. Thus, WDEC is implicitly cal- 
culating a delay based upon 



Figure [i] s hows a diffe rent initial oxygen profile which is 
computed in |Salaris et al. (1997) for a 0.61 Mq white dwarf 
model. This profile was obtained by considering nuclear re- 
action processes in the white dwarf progenitor. Here we use 
a modified algorithm for mixing which reduces to the "com- 
plete mixing" algorithm when applied to an initially flat dis- 
tribution. When a shell crystallizes, we check to see if the 
enhanced carbon content of the innermost fluid shell now has 
more carbon than the shell overlying it. If it does, then we mix 
the two shells and perform the same comparison with the next 
shell farther out, mixing all three shells if necessary. In this 
way, we move outward through the fluid until further mixing 
no longer decreases the carbon content of the fluid between this 
point and the crystallization boundary. This is physically rea- 
sonable, since carbon is, in this sense, "lighter" than oxygen, so 
these layers should be mixed by a convective instability. 

For com pletene ss, we use a third profile taken from Wood] 
(1990) and Wood] (1995). It is designed to be repre sentative of 
C/O profiles calculated in |Mazzitelli & D'Antona (1986) and 
D'Antona & Mazzitelli (1989), who also consider nuclear reac- 



1996; tion rates. Algebraically, it is given by 



X m — 




1.875(^-0.5) 



0.0<4<0.5 
0.5<g<0.9 
0.9 <4< 1.0 



(7) 



(6) 



where q = M^/M* and X c - 1 -X ox . 

Our treatment of the mixing process provides an upper bound 
for the efficiency of this process. If we were to perform a 
more self-consistent calculation, we would compute the Brunt- 
Vaisala frequency for a given chemical composition profile in 
the model and mix those layers which were convectively un- 
stable and whose computed timescales for mixing were shorter 
than the individual timesteps in our evolutionary calculations. 
A n analytical approach t o this more de tailed problem is given 
in tsern et al. (1997) and Mochkovitch (1983). Here we merely 
note that a typical value of \N Z | for a Rayleigh-Taylor unstable 
region in the cores of our models is ~ 10~ 4 , yielding a timescale 
for the mixing instability of ~ 10 2 s, which is clearly shorter 
than the relevant timescales for evolution. 



instead of the expression in equation Q We can remedy this 
situation by an appropriate rescaling of e ps . If we rescale A£,, 
and hence e ps , by 2L,/(L,+L,_i), then equation |6| is transformed 
into equation 0, and we recover the correct age delay due to 
crystallization when implemented in the evolution code. In the 
limit that our timesteps are very small, the above prescription is 
not necessary, but such small timesteps would be computation- 
ally inconvenient, both from a cpu-time standpoint and from a 
numerical convergence standpoint. 

3.3. Consistency Checks 

We use three different initial C/O profiles in our analysis. In 
Figure || we show the oxygen composition in the core both be- 
fore (dotted line) and after (solid line) crystallization has taken 
place. We have taken a homogeneous 50:50 C/O initial distribu- 
tion and assumed complete mixing of the overlying fluid layers 
as crystallization takes place. This should place an upper limit 
on the effect which phase separation can have on any particular 
model. We note t hat the composition pr ofile after crystalliza- 
tion assuming the Segretain & Chabrier] (199 3) phase diagram 
agrees well with that given in Chabrier et al. (1993). 



4. A SIMPLE TEST PROBLEM 

As a check of the standard approach to treating phase sep- 
aration, we perfo rmed a simplified treatment of that given in 
Xu & Van Horn (1992), in which they calculate the change in 
binding energy of a zero-temperature C/Fe white dwarf which 
undergoes phase separation. In order to do this, we have writ- 
ten a separate code which implements the equations for a zero- 
temperature degenerate electron gas. Our approach is simpler 
in that we do not include Coulomb effects in our EOS calcula- 
t ions, so our app roach is essentially pure Chandrasekhar theory 
( Chandrasekhai 1939). We do include, however, relativistic ef- 
fects, which Xu & Van Horn are unable to treat. In testing this 
approach to phase separation, we compute the energy released 
due to phase separation in two different ways. First, we directly 
compute the global change in the binding energy. Second, we 
use the expression for the local energy release and integrate this 
over the mass of the model, as given in equation ||. To fur- 
ther simplify things, we have taken the initial state to be one in 
which the distribution of Fe and C is uniform throughout the 
model, and we have taken the final state to be a pure Fe core 
surrounded by a pure C mantle. 
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FIG. 2. — A comparison of the energy released during crystallization from a static calculation (line) with that from a self-consistent evolutionary calculation (filled 
dots). The error for the total energy released at complete crystallization is less than 0.5%. 
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ture, as shown by the solid line (a). After crystallization is complete, the oxygen 
(1993) is used, and by the dashed line (c) if the phase diagram of Ichimaru et al. 
(1988) is used. We have assumed complete mixing of the remaining fluid layers at each stage of crystallization. It is the redistribution of matter from the initial to 
the final profile which results in a net decrease in the overall binding energy of the configuration. This model is for a 0.6 Mq white dwarf. 



FIG. 3. — Before crystallization has occurred, we have asst 
profile is given by the dotted line (b) if the phase diagram of Segretain & Chabriei 
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FIG. 4. — The same as Figure [?J except that flip initial r/n profile (a) is that computed by Salads et al. (1997) for a 0.61 Mq white dwarf model. Curves (b) and 
(c) are the final profiles assuming the Segretain & Chabriei (1993) and Ichimaru et al. (1988) phase diagrams, respectively. Note that the oxygen mass-fraction at 
the very center increases by only about 15% during crystallization in this case, as compared with a 40% increase for the central value in Figure H. Thus, less energy 
is liberated. 



Figure |5J shows the results for differing initial fractions of 
C and Fe. For instance, for an initial 50:50 C/Fe distribution 
we calculate an energy release of about 1.9 • 10 48 ergs, with 
less than a 5% relative error between the two methods. Even 
this small amount of error decreases as we approach a pure 
Fe or C initial state. This is because the density and composi- 
tion changes before and after phase separation are now smaller, 
which makes the local calculation more accurate. For instance, 
if the model is 99% C uniformly distributed initially, then af- 
ter phase separation most of it (99% in fact) is pure C. The 
contrast between 99% and 100% is small enough that the lo- 
cal density and composition changes are also small (SXq <J 1%, 

^ <; 0. 1 %), which means that the approximation involved in 
making the infinitesimal variations in equation |^ finite is more 
accurate. We note that it is possible to perform such a simplified 
treatment for a C/Fe white dwarf model and still obtain mean- 
ingful results, while for a C/O model it would not be possible. 
This is because \i t = 2.0 for both C and O, while fi & = 2.15 for 
Fe. Thus, ignoring Coulomb effects, C and O have identical 
equations of state, while C and Fe are still nontrivially different 
in this approximation. 

The results of this test problem (Figure ||) convince us that by 
applying equation |[ we are correctly calculating the change in 
the binding energy of the configuration, and thus the amount 
of thermal energy which has been liberated from structural 
sources. This shows us that the overall approach to this prob- 
lem which we use here, and which has been used in the past, is 
sound and accurately describes the physics of phase separation. 

5. RESULTS 

5.1. 0.6 Mq White Dwarf Models 



In Table [j] we give an evolutionary listing of our fiducial se- 
quence (other sequences are available from the author upon re- 
quest). This sequence is more than just a convenient reference 
model for the rest of our calculations. Given the observed peak 
of the masses of isolated whit e dwarfs in the vicinit y of 0.6 M© 
( Weidemann & Koesterj 1983;|Weidemann & Yuan 1989; B erg- 



eron et al. 1995; ^amontagne et alj 1997), this model will be 
the most useful in our comparisons with the white dwarf popu- 
lation as a whole. For the surface layer masses , we have taken 
Mne/M+ = 10" 2 and M H /M* = 10" 4 as in |vVood| (1995). We ex- 
plore the effect of different surface layer masses later in this 
section. 

We now wish to consider the effect of phase separation on 
actual evolutionary sequences. We compute this effect in two 
different ways. Using the first method, we compare sequences 
in which the physics of phase separation is included with those 
in which it is not. Taking two such sequences, we first per- 
form a spline fit for each sequence's age over a fixed luminos- 
ity grid, and then we calculate the difference in ages at each 
luminosity. The results are given by the solid line in Figure ^ 
and indicate an age delay at complete crystallization of about 
1.5 Gyr; this is for an initially homogeneous 50:50 C/O pro- 
file with an assumed metallicity of Z = 0.000 in the opacities. 
We note that all the age delays computed here are for complete 
crystallization of our models, and hence represent upper limits 
to the possible age delays. Operationally, our models are almost 
completely crystallized near the observed luminosity turndown 
at log L^/Lq ~ -4.5, so that there is at most a 1% change in our 
calculated age delays if we only consider models which have 
not yet cooled beyond this point. 

The other method involves applying equation |] to a sequence 
undergoing phase separation, which is shown by the dotted line 
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X c (Carbon Mass Fraction) 

FIG. 5. — The upper panel shows the energy released due to phase separation as a function of X c , the carbon mass-fraction. The model is a C/Fe mixture computed 
assuming pure Chandraaekhar theory. The solid line comes from a direct calculation of the change in binding energy, and the dotted line is obtained from the 
application of equation H The lower panel shows the percent error between these two methods The total mass of the model is set tn f) b^Afi M c as in Xn & Van 
Horn (1092). ' 



Table 1 

Fiducial White Dwarf Cooling Sequence without Phase Separation 



log L/L Q 


log Age (yr) 


logT e 


logoff 


logfl+ 


log L„/L Q 


M x tal/M* 


1.0000 


5.917 


7.950 


4.869 


9.131 


1.298 


0.000 


0.6000 


6.158 


7.891 


4.793 


9.084 


0.946 


0.000 


0.2000 


6.376 


7.843 


4.710 


9.049 


0.581 


0.000 


-0.2000 


6.592 


7.798 


4.623 


9.024 


0.193 


0.000 


-0.6000 


6.855 


7.737 


4.533 


9.004 


-0.230 


0.000 


-1.0000 


7.204 


7.660 


4.440 


8.989 


-0.886 


0.000 


-1.2000 


7.429 


7.604 


4.394 


8.982 


-1.356 


0.000 


-1.4000 


7.674 


7.531 


4.346 


8.976 


-1.956 


0.000 


-1.6000 


7.903 


7.448 


4.299 


8.971 


-2.670 


0.000 


-1.8000 


8.097 


7.360 


4.251 


8.967 


-3.449 


0.000 


-2.0000 


8.264 


7.273 


4.203 


8.963 


-4.248 


0.000 


-2.2000 


8.413 


7.187 


4.155 


8.959 


-5.042 


0.000 


-2.4000 


8.550 


7.103 


4.106 


8.956 


-5.804 


0.000 


-2.6000 


8.679 


7.021 


4.058 


8.953 


-7.279 


0.000 


-2.8000 


8.805 


6.940 


4.009 


8.951 


< -10.00 


0.000 


-3.0000 


8.930 


6.860 


3.961 


8.948 


< -10.00 


0.000 


-3.2000 


9.055 


6.778 


3.912 


8.946 


< -10.00 


0.000 


-3.4000 


9.182 


6.694 


3.863 


8.944 


< -10.00 


0.000 


-3.6000 


9.317 


6.613 


3.814 


8.942 


< -10.00 


0.059 


-3.8000 


9.497 


6.518 


3.765 


8.939 


< -10.00 


0.379 


-4.0000 


9.704 


6.350 


3.717 


8.935 


< -10.00 


0.804 


-4.2000 


9.825 


6.176 


3.668 


8.933 


< -10.00 


0.961 


-4.4000 


9.900 


6.022 


3.618 


8.933 


< -10.00 


0.989 


-4.6000 


9.953 


5.891 


3.569 


8.932 


< -10.00 


0.990 
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FIG. 6. — The solid line is a self-consistent calculation of the age difference between two 0.6 Mq white dwarf evolutionary sequences with Z = 0.0, one of which 
is undergoing phase separation. The dotted line is the result of applying equation M to the evolutionary sequence undergoing phase separation, which yields an 
asymptotic value for the age delay of ~ 1.4 Gyr. At complete crystallization (logL/tQ ^l.— 4.6), the value given by the direct evolutionary calculation is within 5% 
of this, indicating that the basic physics which is operating is well-described by equation hi 



in Figure This yields an asymptotic value for the age delay 
of 1.38 Gyr, which is within 5% of the age difference com- 
puted with the first method. This result indicates that the basic 
physics which is operating is well-described by equation |], i.e., 
the energy being released by phase separation is mostly being 
radiated in a given timestep. For the remainder of the results 
quoted here, the age delays have been calculated using this sec- 
ond method (equation ||), since this proves to be more accurate 
for the cases involving smaller energy releases and age delays. 

We now study the effect of the initial composition profile on 
the age delays. We use three different profiles: one which is a 
hom ogeneous 50:50 mix (Figure ^|), one calculated by S alaris 
et al. (1997) (Figure and one given by equation [7[ Our re- 
sults are summarized in Table |[ where the columns labeled SC 
and IIO indicat e that we have used the phase diagrams of S egre- 
tain & Chabrier (1993) and Ichimaru et al. (1988), respectively. 
Near the centers of these models, we found that the initial/final 
oxygen mass-fraction changed by only about 15% in the ini- 
tially stratified case in Figure ||, as compared to 40% in the 
homogeneous case in Figure |[ Because less matter is being re- 
distributed in the initially stratified case, we would expect less 
energy to be released as a result. Using the phase diagram of 
Segretain & Chabrier (1993) applied to a 0.6 Mq white dwarf 



model, we find that in the homogeneous case 2.38 • 10 ergs are 



released whereas in the initially stratified case in Figure [4| only 
1 .03 ■ 10 46 ergs are released. These energies result in age delays 
of 1.38 Gyr and 0.62 Gyr, respectively. Thus, the initial com- 
position profile has a large effect on the calculated age delays. 
In addition, the Ichimaru et al. (1988) phase diagram produces 
smalle r composition changes and hence smaller values, reduc- 
ing the Segretain & Chabriei age delays by approximately one- 
third. 



We now consider the effect of a nonzero metallicity in the 
opacity tables. The effect of varying the metallicity from 
Z = 0.000 to Z = 0.001 results in a change of less than 0.016% 
in the energies released, and is barely detectable numerically. 
The main effect of changing the metallicity is to affect the lu- 
minosity range at which the phase separation energy is released, 
which in turn affects the age delay, fdeiay For both the homo- 
geneous and stratified case, the average luminosity during crys- 
tallization changes by less than 3% as Z is varied from 0.000 to 
0.001, and hence fdeiay also changes by less than 3%. Thus, the 
age delay is essentially insensitive to the metallicity assumed 
for the opacities. 

Finally, we summarize the effect of different surface layer 
masses in Table (| For M He /M+ = 10" 3 and M H /M+ = 10" 5 
(compl), we find maximum age delays of 1.45 Gyr, and for 
Mfte/M-t = 10" 4 andMn/M* = 10" 6 (comp2), our maximum cal- 
culated age delay is 1 .56 Gyr. These values represent increases 
of 5% and 13%, respectively, over the age delays calculated in 
our fiducial model. For clarity, we note that these calculations 
are for the age differences introduced by phase separation alone 
at these new surface layer masses; the white dwarf ages them- 
selves change significantly with He layer mass, which produces 
a decrease in the calculated ages (without including phase sep- 
aration) of ~ 0.75 Gyr for each order of magnitude increase in 
Mue- Again, we find that varying the metallicity in the opacities 
has a small effect on these numbers, at only the 1% level. 

5.2. The Mass Dependence 

The mass of the white dwarf model affects the process of 
phase separation in two main ways, as is illustrated in Figure 
Uj First, a more massive white dwarf has a higher gravity so 
that more energy is released by the subsequent rearrangement 
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Table 2 

Age Delays for O.6M Models 



SC 



Delay (Gyr) 



IIO 



1.38 
0.62 
0.30 



0.99 
0.39 
0.20 



Initial Profile 



50:50 ho mogeneous 



stratified ( Salaris et al 
stratified ( |Woodj 1995) 



1997) 



Table 3 

Age Delays for 0.6M q Models with Different Surface Layer Masses 



Initial Profile 



compl 



Delay (Gyr) 



comp2 



SC 


IIO 


SC 


IIO 


1.45 


1.04 


1.56 


1.12 


0.66 


0.42 


0.71 


0.44 


0.32 


0.21 


0.34 


0.23 



50:50 ho mogeneous 



stratified (Salaris et al 
stratified ( |Wood| 1995) 



1997) 



of matter. Second, the luminosity at which crystallization oc- 
curs is higher for a more massive white dwarf, which tends to 
lessen the age-delay for a given energy release. For example, 
even though the total energy released in a 1.2 M Q model in- 
creases by a factor of ~ 10, the average luminosity increases by 
a factor of ~ 30, and hence there is a net decrease in the time 
delay relative to the 0.6 M Q sequence. 

The competition of these two effects suggests that there may 
be a mass for which there is a maximum age-delay, for a fixed 
composition profile. This is indeed the case, as is demonstrated 
in Figure ^[ We find that the 0.6 M Q white dwarf models have 
the maximum age- delays for a giv en composition profile (this 
was also found by Segretain et al. 1994). The calculated age 
delay is only weakly dependent upon the metallicity, as can be 
seen from the small difference between the solid and dashed 
curves. It is strongly dependent upon the initial profile, how- 
ever, which can decrease the energy release, and hence the age- 
delays, by a factor of three or more, as is shown in Figure ||. 

From the preceding calculations we find that the two most 
important factors influencing the magnitude of the age delays 
introduced by the physics of phase separation are the mass of 
the white dwarf model and its initial C/O profile. Because the 
mass range o f observed white dw arfs is strongly peaked around 
0.6M Q (e.g., I^amontagne et al. 1997), we find that the age de- 
lay we calculate is near the maximum possible with respect to 
this parameter. In terms of the initial C/O profile, however, the 
situation is revers ed. For a 0.6 1 M B white dwarf model, the pro- 
file calculated by Salaris et al. (1997) reduces the age delay by 
a factor o f ~ 2 from the 50:50 homogeneous case. Using the 
profile of |Wood| ( 1995), which is based on results from Mazz- 
itelli & D'Antona (1986) and p'Antona & Mazzitelli| (1989), 
the reduction factor is ~ 5. 



Gyr. 



If we take as our best guess the initial profile of |Salaris et al 
(1997), assume a O.6M white dwarf model with Mu e ]M± = 
10" 2 and M H /M t = 10" 4 , and use the |Segretain & Chabriei 



6. CONCLUSIONS 



We find a maximum age delay of ~ 1 .5 Gyr due to phase sep- 
aration for our fiducial white dwarf model (M* = .6M<^) and a 
best guess age delay of ~ 0.6 Gyr. Salaris et al. (1997) have 
recently calcu lated a value of ~ 1 Gyr, using the evolution- 



ary models of Wood & Winget (1989). If we scale their value 
to our present models (assuming an average luminosity during 
crystallization for their models of logL/L Q ~ -4.1), then we 
obtain 0.75 Gyr, which is in basic agreement with our estimate 
of 0.62 Gyr. The differences in these models are mainly due 
to the different surface layer masses adopted; more recent as- 
teroseismological analyses of the class of DA's suggests that 
the appropriate surface la yer masses are Mu e /M+ ~ 10~ 2 and 
M H /M+ ~ 10" 4 (|Clemens| 1993, 1995), and these are the values 



(1993) phase diagram, then we obtain an age delay of ~ 0.6 



which we have assumed. 

The most important factors influencing the size of the calcu- 
lated age delay are the total stellar mass and the initial composi- 
tion profile. We find the largest age delays occur in models with 
masses of ~ 0.6M Q , near the peak in the observed white dwarf 
mass distribution. The best current theoretical initial C/O pro- 
file produces models with smaller age delays, of ~ 0.6 Gyr. In 
addition, if w e use the phase diagra m of Ichimaru et al. (1988) 
instead of the Segretain & Chabrier (1993) phase diagram, then 
our age delays are reduced by about one-third. We note that 
the prescription which we have adopted for the mixing during 
crystallization provides an upper bound for the efficiency of this 
process, and hence a maximum for the age delay. More realistic 
treatments of the mixing process may reduce the age delay. We 
find that varying the opacities (via the metallicity) and varying 
the surface layer masses has only a small effect (<; 10%) on the 
calculated age delays. 

Our calculations do not take into account the possible age 
delays introduced by the phase separation of heavier trace- 
element species such as 22 Ne, which may produce significant 
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FIG. 7. — Phase separation energy and average luminosity as a function of mass, using the Segrer ain & rhahrip r phase diagram. The solid curves are for initially 
homogeneous 50:50 C/O mixture s the dashed curves are for the stratified C/O profile of Figure M ( jSalaris et al.| 1997), and the dotted curves are for the stratified 
C/O profile of equation [7 ( Wooc 1995). The upper panel shows the total phase separation energy released as a function of total stellar mass, and the lower panel 
shows the average luminosity during the crystallization process, also as a function of total stellar mass. 
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FIG. 8. — Age delay due to phase separation during crystallization as a function of total mass of the white dwarf model. Curve (a) corresponds to a 50:50 
homogeneous initial C/O profile, while curves (b) and (c) are the initial profiles shown by the solid lines in Figures W and M, respectively. The solid lines are for 
zero metallicity opacities and the dashed lines are for Z = 0.001, which shows that our result has little metallicity dependence. All models have Mne/M* = 10~ 2 and 
M H /M* = 10-*. 
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age delays of 2-3 Gyr ( Segretain et al. 1994; frlernanz et al 
1994). These species would arise from the initial abundance 
of metals in the main-sequence stars which later evolved into 
white dwarfs. This effect may only be important for Population 
I stars, however, and would not therefore affect the calculated 
ages of the cool white dwarfs which populate the turndown in 
the WDLF, since these w hite dwarfs wer e formed very early in 
the history of the Galaxy (Hernanz et al. 1994). 

In the context of Galactic evolution, age estimates for the 
olde st Galactic globular clusters range from 13-16 Gyr (P ont 
et al. 1998) to 11.5 ~ 1.3 Gyr ( piaboyeret al.| 1998), and de- 
pend on a variety of parameters. In addition, a 4 to 6 Gyr delay 
is expected between the formation of the globular clusters and 
that of the Galactic thin disk (e.g., Burkert, Truran, & Henslei 



1992; phiappini, Matteucci, & Grattorj 1997), while the ob- 
served white dwarf luminosity f unction gives an age estimate 
for the thin disk of 9.5to \ Gyr (|Oswalt et ak| 1996), without 



including the effect of phase separation. Using the above num- 
bers, we see that phase separation could add anywhere from 
to 3 Gyr to the white dwarf ages and still be consistent with the 
overall picture of Galaxy formation. Our calculated maximum 
value of <; 1.5 Gyr fits within these bounds, as does our best 
guess value of ~ 0.6 Gyr. 
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